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We study the even-parity I — 2 perturbations of a Schwarzschild black hole to second order. 
The Einstein equations can be reduced to a single linear wave equation with a potential and a 
source term. The source term is quadratic in terms of the first order perturbations. This provides 
a formalism to address the validity of many first order calculations of interest in astrophysics. 



\& ' Black hole perturbation theory has been a ubiquitous tool in the analysis of astrophysical situations without 
0^ i symmetries. It has played an important role in the study of the expected gravitational radiation from processes like 
the infall of matter into a hole jij , and more recently the collision of two black holes Q . For black holes, perturbations 
to linear order in a dimensionless expansion parameter e can be described through a simple scalar wave equation in 
two spacetime dimensions, and yet can describe complex situations of no particular symmetry. In contrast, a full 
numerical simulation of the nonlinear Einstein equations is at present prohibitive for most situations of physical 
'™~ 5 ■ interest. 

In spite of these advantages there is an important limitation in the use of linearized theory: there is no information 
within linearized theory to determine its range of applicability, i.e., to determine what values of e are sufficiently small. 
' Largely due to this last problem, linearized theory has usually been limited in use to situations departing only very 
slightly from an undisturbed black hole. However, recent comparisons between linearized theory and full numerical 
' simulations for the head-on collisions of two black holes have shown that linearized theory can work remarkably 
well in domains in which it would be supposed to fail. The qualitative explanation for this is that strong nonlincaritics 
■ near the horizon can be absorbed by the black hole and therefore not affect the outgoing radiation. 

In order to promote linearized theory to a serious tool for predicting astrophysical answers and in particular to 
provide benchmarks for difficult full numerical simulations, we need a reliable measure of the errors in a perturbation 
result for a given value of e. One can try to construct simply implemented analytic a priori indices of the validity of 
perturbation theory, e.g., the violation of the exact Hamiltonian constraint by the linearized data or the violation 
of the linearized Hamiltonian constraint by the exact initial data [Q . Such measures are interesting because of their 
simplicity, but are not unique and worse, they are not guaranteed to work in all situations. For instance, their success 
£~ I ■ for head-on collisions of two black holes gives no reliable information about their value for less symmetric collisions. 

There is, in fact, only one generally reliable index of the accuracy of linearized perturbation theory: to calculate 
answers to the next order in e. Where the higher-order results and the linearized results differ by (say) 10% is a point 
at which one has some confidence that either answer is accurate to within around 10%. Higher order perturbation 
calculations seem at the outset to be an obvious extension of the familiar techniques of linearized perturbation theory, 
with a guarantee of simple results. It turns out — and this is one of the major points we wish to make — that there 
are some subtle issues of gauge that must be recognized. Despite these issues, and despite some necessarily lengthy 
calculations, second order computations are far more easily done than the development of full numerical relativity 
codes. And the difficulty of the higher order calculations is a necessary price well worth paying for the advantages 
provided by reliable perturbation theory results. 

The basic structure of higher order perturbation theory starts with the same basic formalism as linearized theory. 
We write the metric as g$ + eg$ + e 2 gffi + C(e 3 ) where g^u, the "background," is a known solution of the vacuum 
field equations. We write the sourceless Einstein equations as 

G(9$+e 9 $ + <?9$) = 0, (1) 

where Q represents the actions of taking partial derivatives and algebraic combinations to form the components of the 
Einstein tensor. If we expand ([!]) in e, the term of order e° automatically vanishes if g$ is a background solution. 
The terms first order in e can be written in the form 

eAff#) = 0, (2) 

where £ is a set of differentiations and combinations with details that depends on gffl . These operations are all linear 
in g^J , and eqs. (0) constitute linearized perturbation theory, a reduction of the problem to a system of linearized 
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equations. If the background is the Schwarzschild geometry, it is well known that these linear equations can be 
decomposed into multipole moments, can be separated into two sets of independent functions (the even-parity and 
odd-parity perturbations), and can be rearranged into two simple linear wave equations in the variables r,t: the 
Regge- Wheeler || equation for the odd-parity perturbations, and the Zerilli || equation for even-parity. 

The part of ([!]) that is proportional to e 2 has two kinds of terms. There are terms that are linear in gffi, and terms 
that are quadratic in g$ . It is clear that the former terms occur in precisely the same form as do the g^J terms in 
(|J). The set of e 2 terms can then be written as 

e 2 £(g$)=e 2 J(g$,g$), (3) 

where J is quadratic in the first order perturbations. In solving for the second order perturbations, one treats the 
first order perturbations as already known, so J plays the role of a source term in (|^). 

The C operator in (||) is precisely the same operator as in (||), so for each linearized theory equation for g$ there 

(2) 

is a corresponding equation for g^J , differing only in the presence of a "source" term. The very same manipulations 
that lead to the Regge- Wheeler and Zerilli equations, must therefore lead to Regge- Wheeler and Zerilli equations of 
precisely the same form for second order perturbations, except for the presence of source terms. We therefore have a 
guarantee at the outset that, as in first order theory, we can derive a simple wave equation for the perturbations. 

The form of the metric functions, to any order in e, depend on how the spacetime coordinates are defined, to that 
order in e. We assume that we have a set of coordinates t, r, 9, 4> in terms of which the metric, to zero order in e, takes 
on the standard Schwarzschild form. The background metric gffl then is the Schwarzschild metric. We can transform 
to new coordinates x a with e-dependent coordinate transformations 
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al - x a + eXfaixf) + e 2 Xf 2) { X P) + 0(e A ) . (4) 



The form of the background metric gffl is unchanged by the transformation, and we call such a coordinate change 



a gauge transformation. If the first order perturbation functions g^v are changed by the coordinate change (i.e., if 
{x 13 ) 7^ 0) we call it a first order gauge transformation. Note that a first order gauge transformation will in general 

(2) 

also change the form of the second order perturbations g^J . We can also make gauge transformations which leave 
g$ invariant and which change gffi (i.e., transformations with X^{x^) — 0, but for which X? 2 Ax^) ^ 0). These we 
call second order gauge transformations. 

We will want to take advantage of the freedom to choose coordinates in order to impose some simplifying special 
conditions on the metric perturbations. As an example, for even parity perturbations we will at one stage be setting 
g r e = through second order in e. To do this one makes a first order gauge transformation to set g^J = 0. This 
transformation will affect the second order perturbations in some way that depends on the details of the coordinate 
transformation. One next makes a second order transformation to set g^ 2 g = 0. This leaves g^ g = 0, since a second 
order transformation cannot change the first order perturbations. 

In higher order perturbation calculations new questions can arise, such as what it means for a quantity to be gauge 
invariant. It turns out to be useful always to think of coordinate transformation as a sequence of distinct steps, as in 
the example above. In this way we see that a quantity can be "second order gauge invariant" (invariant only under a 
second order transformation) , "first order invariant" (invariant only to first order under a first order transformation) 
or "first and second order invariant" (invariant up to second order for a general transformation (Q)). 

In this paper we take the first important steps in using second order perturbation theory as a tool; we provide the 
formalism for calculating the second order contribution to outgoing radiation. We will, however, make important, and 
yet physically sensible restrictions in our formulation. First, we will restrict attention to the even-parity t = 2, to = 
second order contributions. This is justified by the fact that most radiation processes are dominated by quadrupolar 
radiation, and we are primarily interested in the "error measure" on the first order quadrupole calculation. We give 
here the even-parity formalism, but the odd-parity equivalent follows the same pattern and is in fact simpler. The 
restriction to to = is for simplicity only, the generalization to to ^ is immediate. Our second, and more subtle, 
restriction is on the contributions to the "source" term. Since the source is quadratic in the first order perturbations, 
there is an infinite number of first order multipoles contributing to the £ — 2 projection of the source. We restrict 
attention only to the first order 1 — 1 contributions. There are several reasons for this: (i) It is primarily a simplifying 
assumption; the £ = 2 projection of the complete source is straightforward to compute either numerically or as an 
infinite series, but would introduce unnecessarily distracting complications here, (ii) As a practical matter, it is 
plausible that in most practical situations the source will be dominated by the first order quadrupole terms, (iii) The 
specific case to which we will first apply this formalism (details to be published elsewhere) is the head-on collision 
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from initially small separation in which the initial separation is the expansion parameter. In that case it can be 
shown that the first order perturbations have only an I — 2 contribution. 

In terms of the usual Schwarzschild-like coordinates t, r, 9, (f), and following the Regge- Wheeler ^ prescription and 
notation, we write the general form of the £ — 2, even parity, perturbation of the spherically symmetric black hole 
metric in the form, 
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where P 2 (9) = 3(cos 2 9 - l)/2, P!,(9) = dP 2 (9)/d9, and P%(9) = d 2 P 2 (9)/d9 2 , and where the functions H, h, K and G 
depend on t and r only. Just as in J5|, which was restricted to linearized theory, we may impose the Regge- Wheeler 

gauge conditions through second order, = = G^ = 0, i = 1, 2, by the two-step process described above. One 
can show, (details will be given in a lengthier paper), that given an arbitrary perturbation of the form (|^), there exists 
always (locally) a uniquely defined gauge transformation that takes the metric to the Regge- Wheeler gauge. (One 
can, in fact, write these Regge- Wheeler functions, both to first and second order, explicitly in terms of perturbations 
in an arbitrary gauge, and through these expressions view the Regge- Wheeler perturbations as gauge invariant.) 

With the linearized equations in the Regge- Wheeler gauge, Zerilli || assumes time dependence e %UJt (i.e., makes 
a fourier transform), and works with functions of u) and r, but the process of deriving a single wave equation can 
be done equally well in terms of the original r, t variables. One can repeat this Zerilli process, step by step, with 
the second order equations, in the Regge- Wheeler gauge, in which the only new feature is the inclusion of quadratic 
"source" terms. One finds a total of seven nontrivial Einstein equations for the four second order Regge- Wheeler 
perturbation functions Hq 2 \ h[ 2 \ K^ 2 \ These equations are linear in the second order functions, but quadratic 
in the first order perturbations. One of these takes the form 
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and may be used to eliminate H^p from the other equations. (We also have = from the first order 

equations). We are therefore left with six equations for the three remaining second order functions, but it turns 
out that we may choose three of these equations so that they can be obtained by compatibility of the other three. 
Moreover, from one of Einstein's equations we find, 



H^\ t = ~K^\ t +[(1 - 2M/r)H[ 2 \ r +(1/14)[2(^«) 2 + 3(^ 1} ) 2 - 3(i^ 1} ) 2 



(7) 
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(a comma indicates partial derivative), which may be used to eliminate Hq (or rather its time derivative) in the 
remaining two independent equations. Following Zerilli we introduce now a pair of functions x(^ r ) an d 7t(t,r), 
related to K^(t,r), and H (2 \t,r) by 



(2) _ 6(r 2 + Mr + M 2 ) 
'* r 2 (2r + 3M) 
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By substitution of these relations in Einstein's equations we find, 
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* = (l - *f ) X,r 4^ (^Wrf -2^)^, -rtf^ ^ 



while the function x(t, r) satisfies the single second order differential equation 



<9i 2 9r* 



+ V(r) X (t,r) =5 RW . (11) 



where r* = r + 2Mlog(r/2Af — 1) is the "tortoise" radial coordinate, and where V(r) is the £ — 2 Zerilli potential, 
identical to that for the linearized theory (see || ) . The "source term" <Skw is a quadratic expression in the first order 
perturbations. 

As guaranteed by the procedure, we have arrived at a Zerilli equation for a second order Zerilli function, differing 
from the Zerilli equation of linearized theory only in the presence of the source term <Srw quadratic in the first order 
perturbations. (The explicit expression for S^iw is somewhat lengthy and will be given elsewhere.) The second order 
Zerilli equation ( jll| ) is then the core of second order perturbation theory. We have verified, as a consistency check, 
that these functions satisfy the second order Einstein equations, (or rather, their time derivatives 0), provided that 
the first order Einstein equations and (|ll| ) are satisfied. In principle, one first solves the first order Zerilli equation, 
finds the first order perturbations, constructs the source term <Srw, and solves (|i~l| ) for \. In practice, however, 
there are complications. These are evident in the fact that the source term in ([L 1|) diverges at large radius. This 
divergence does not indicate singular physical behavior; the second order radiation could in principle be extracted 
from x computed in this way, but the procedure would, at best, be computationally inefficient. In connection with 
this, two issues concerning gauge choice should, in particular, be noted. First, it must be realized that most of the 
development of the second order Zerilli equation (|ll) did not require that the gauge choice be Regge- Wheeler for the 
first order perturbations; only second order Regge- Wheeler restrictions are really needed. If we had made another 
gauge choice, or no gauge restrictions at all, only the quadratic terms in (Q), ([n]), and (|ll|) would change. The "second 
order Zerilli equation with source" is, therefore, not unique. Another viewpoint on this is that we could start with 
( |j~T| ) and introduce a new Zerilli function through 

X = Xncw + quad (12) 

where quad is any expression quadratic in the first order perturbations. The new Zerilli function Xncw would then sat- 
isfy ([ll]), but with a modified source term. The choice of the Regge- Wheeler first order gauge is made for convenience. 
Expressions in this gauge usually turn out to be most compact. 

The second gauge issue of note concerns the computation of radiation. To first order all information about radiation 
is carried in a "Zerilli function" ip which we define, in the Regge- Wheeler gauge, with the notation introduced in (j|), 
by 

r(r-2M) ( m dK^\ r „ m 

^ = oL , zrk N 3 + -k {1) ■ (13) 



3(2r + 3M) V dr 

Several different normalizations for the "Zerilli function" have been used in the literature, so it is important to specify 
the relationship of our choice to others. If there are no (first order) sources, the definition in ([l3]) can be shown, for 
£ = 2, to have the same appearance as the definition of the wave funtion Rlm of Zerilli. That is, ( |l3|) specifies the 
same combination of perturbation functions and so forth. But in our definition (|s|) of the metric perturbation 
functions, we have expanded in Legendre polynomials, whereas Zerilli does his multipole expansions in terms of 
spherical harmonics. As a result the actual relationship between our ip and Zerilli's R is 



(i4) 

Another normalization that is of importance is the normalization used by Cunningham et al. [^) , which agrees with 
that of Zerilli. The normalization used in and [Eof , called import in the latter reference, is related to ip used here by: 



f2£TT n (£~2)l , 
V 4tt (£ + 2)r port ' 



^ = \/^^2)— ^ pert . (15) 



4 



The computation of graviational radiation power is done in a coordinate system that is asymptotically flat. In this 
system, which we will denote with tildes, radiation information is carried by the perturbations in g#6),<?00 and gee- 
The function G, in the notation of (||), falls off in this gauge as 1/r, to all orders in e, and rG can be thought of as 
the amplitude of the even parity gravitational wave. In terms of our definitions, The first order part is given by 
C?,j = ip(t, r), t jr + 0(l/r 2 ) and the first order radiated power can be shown to be 

Power = — j £ + 2 j| («/>) 2 , (16) 

16 {£-2)1 21 + V V > ' y 1 

where ip is the time derivative of ip. For I = 2 the result is Power = (3 / '10) (etp) 2 . If we want the gravitational wave 
amplitude to second order we must calculate eG^ + e 2 G^ 2 \ 

The computation of & 2 ^ can be approached in several ways. One could, in principle, transform to a gauge which is 
asymptotically flat to first order, so that the source term <Srw is replaced by the appropriately modified source term 
iS ra d. One then solves for Xi from it computes the second order Regge- Wheeler metric perturbations, and then does 
a second order gauge transformation to a second order asymptotically flat gauge. In this gauge one reads off G^ 2 \ 
In practice the same thing can be accomplished more conveniently with a transformation of the form (|l2]). We start 
with ( pT| ) and with a gauge which has Regge- Wheeler restrictions to first and second order, and we introduce a new 
function x(t, r )rcn, a sort of "renormalized x," given by 
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Simple replacement in (|ll[) gives us an equation for x(i,r) ron with the same form as (|l l|) , but with a source term 
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where A = (2r + 3M), where /x = (r — 2M), and where V' is the solution of the first order Zerilli equation. This 
"renormalized" source dies off asymptotically for large r, and we find that Xron behaves asymptotically as a function 
of t — r* only. We will postpone a detailed proof to a lengthier publication. 
Most important, it can be shown that x ren and the radiation are related by 
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+ 0{l/r 2 ) (19) 



and we see that Xren determine completely, (and in a numerically meaningful way), the asymptotic behavior of & 2 \t- 
It seem appropriate, then, to refer to Xron as the second order Zerilli function. 

To conclude, we present the expression of the gravitational radiation power, accurate to second order in e, from 
which the total radiated energy is easily obtained, in terms of the first and second order Zerilli functions. To compute 
it, we just write the Landau-Lifshitz pseudotensor in terms of the metric perturbations in the asymptotically flat 
gauge (see || for details) and substitute the expression for the manifestly asymptotically flat form of the metric 
perturbation (|l9|). The result is 
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The formalism, finally, consists of the following steps: (i) From a solution to the initial value problem, one extracts 
first order perturbations of three geometry and extrinsic curvature, and from them computes Cauchy data for the first 
order Zerilli function ip (ii) The Zerilli equation is then solved numerically for ip in the t,r region of interest, (iii) 
The solution for ip is used to compute iS rcn in ( |l§| ) . (iv) From the initial value solution one next uses the definition of 
X in the Regge- Wheeler gauge 



from (||) and and the other second order Einstein equations, at the initial hypersurface, to find the initial value 
of \, and its time derivative, (v) Next (|17j ) and its time derivative are used to find the initial Xrcn and x ron (vi) 
The renormalized Zerilli equation [eq. (|ll[) with the source iS ron ] is then solved numerically for Xrcn- (vi) Finally, 
the outgoing radiation power is computed from (p0|). We will be pursuing several applications of this formalism in 
subsequent publication. 
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